#include <stdio.h>
#include <stdlib.h>
#include "linsyst.h"
#include "logtime.h"

#define LOGFILE "mjr-cpu.log"

void transformaestrela(syst * lsys);
void jacobi_richardson(syst * lsys);

int main(int argc, char** argv){
	int i, j;
	float teste, el;
	syst * lsys;

	if(argc != 2){
		printf("Usage: %s FILE\n\n", argv[0]);
		return 1;
	}

	initlogtime();

	// Carrega arquivo
	lsys = loadsyst(argv[1]);
	if(lsys == NULL) return 1;

	regtime(LOADFILE);

	// Aplica o Metodo de Jacobi-Richardson
	transformaestrela(lsys);
	jacobi_richardson(lsys);

	regtime(METJACOBI);
	setit(lsys->max_it);


	// Exibi resultado
	printf("Iterations: %d\n", lsys->max_it);
	
	i = lsys->ref;
	el = lsys->el_digref;
	teste = 0;
	for(j = 0; j < lsys->size; j++){
		if(j != i) teste += MAT_A(lsys,i,j)*(lsys->vetX0[j])*el;
		else teste += (lsys->vetX0[j])*el;
	}

	printf("RowTest: %d => [%f] =? %f\n", i, teste, lsys->vetB[i]*el);

	destroysyst(lsys);


	//registra execucao no arquivo de log
	savelog(argv[1], LOGFILE);

	return 0;
}


void transformaestrela(syst * lsys){
	int i, j;
	float el_dig;

	for( i = 0; i < lsys->size; i++){
		el_dig = MAT_A(lsys,i,i);

		 // A -> A*
		MAT_A(lsys,i,i) = 0;
		for( j = 0; j < lsys->size; j++)
			MAT_A(lsys,i,j) /= el_dig; 

		// B -> B*
		lsys->vetB[i] /= el_dig;
	}
}

void jacobi_richardson(syst * lsys){
	int it, i, j;
	float aux, max_dif, max_x;
	int size;

	size = lsys->size;

	for( i = 0; i < size; i++) // X0 = B*
		lsys->vetX0[i] = lsys->vetB[i];

	for( it = 0; it < lsys->max_it; it++){ // iteracao
		for( i = 0; i < lsys->size; i++){  // calculo de Xk+1[i]
			lsys->vetX1[i] = lsys->vetB[i];
			for( j = 0; j < lsys->size; j++) // A*[i][j] x Xk[j]
				lsys->vetX1[i] -=  (lsys->matA[size*i+j])*(lsys->vetX0[j]);
		}
		max_dif = max_x = 0;
		for( i = 0; i < size; i++ ){ // erro max( |Xk+1[i] - Xk[i]| ) / max( |Xk+1[i]| ) 
			aux = lsys->vetX0[i] - lsys->vetX1[i];
			if(aux > max_dif) max_dif = aux;
			else if(-aux > max_dif) max_dif = -aux;
			aux = lsys->vetX0[i] = lsys->vetX1[i];
			if(aux > max_x) max_x = aux;
			else if(-aux > max_x) max_x = -aux;
		}
		
		if (max_x != 0 && (max_dif/max_x) <= lsys->prec){
			it++;
			break;
		}
	}

	lsys->max_it = it;
	
}



